Localization of Wannier functions for entangled energy bands 
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A new method for the localization of crystalline orbitals for entangled energy bands is proposed. It 
is an extension of the Wannier-Boys algorithm [CM. Zicovich- Wilson, R. Dovesi, and V.R. Saunders, 
J. Chem. Phys. 115, 9708 (2001)] which is particularly well-suited for linear combination of atomic 
orbital representations of the Bloch waves. It allows the inclusion of additional bands during the 
optimization of the unitary hybridization matrix used in the multi-band Wannier transformation. 
By a projection technique, the proper chemical character is extracted from the Bloch waves and 
compact localized orbitals are obtained even for entangled bands. The performance of the new 
projective Wannier-Boys localization is demonstrated on the low-lying unoccupied bands of trans- 
polyacetylene, diamond and silicon. 



PACS numbers: 71.15.Ap,71.15.Dx 
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I. INTRODUCTION 

Localization of Wannier functions (WF) has attracted 
great attention of scientists in the recent years. Despite 
the success of describing most of the physical phenom- 
ena in crystals in terms of Bloch waves (BW), Wannier 
functions have obvious advantages. To mention a few: 
WFs provide a chemically intuitive picture of the elec- 
tronic structure in crystals, using localized WFs, physical 
quantities such as effective Born charges and spontaneous 
polarization can be evaluated in a very simple w&yi^l 
and they play a central role in many post Hartree-Fock 
electron correlation methodsA^SiL^ 

Several rigorous schemes for the localization of WFs 
have been proposed so far. They fall into two categories, 
those which are based on the variational principle 
and those which are an extension of the Foster-Boys lo- 
calization criterioniiii^ for periodic systemsiiii^i^i^ All 
these methods can only be applied to isolated band com- 
plexes, i.e. a group of bands which are separated from 
the other bands by an energy gap over the entire Bril- 
louin zone. This restriction appreciably confines the pos- 
sible applications of those methods and only a limited 
number of systems can be treated. In particular, the 
energy bands of the unoccupied Bloch waves usually do 
not exhibit any such band gaps. Thus, none of the lo- 
calization schemes developed so far can be applied rou- 
tinely to selectively generate virtual Wannier functions 
as needed, for instances, in our wave-function-based post 
Hartree-Fock correlation methods for valence and con- 
duction bands. 17 ' 18 

Recently, Souza et al. 19 have extended the origi- 
nal Marzari-Vanderbilt localization scheme^ to systems 
with entangled bands. The method is based on a pre- 
selection of optimal Bloch waves having maximal simi- 
larity at neighboring fc-points by minimizing a suitable 
functional. Like the original scheme, the new method is 
especially designed for a plane wave representation of the 
BWs and heavily relies on numerical fc-spacc differentia- 
tions. 

In this paper, we propose a similar extension for the 



Wannier-Boys localization scheme developed by Zicovich- 
Wilson, Dovesi and Saunders»i& This scheme differs in 
many aspects from the algorithm proposed by Marzari et 
alH In particular, it is much better suited for BWs given 
in a linear combination of atomic-like orbital (LCAO) 
representation as employed in many widely used pro- 
gram packages for periodic systems such as Crystal jSS 
Gaussian, 21 DMol 3 )S Nfp-Lmto, 23 or Band24 (the 
periodic variant of Adb 2 ^). 

In Sec. [HI the details of our algorithm are described. 
Then, in Sec. lIIII the results from the localization are pre- 
sented for irans-polyacetylene (tPA), diamond and sili- 
con, and some concluding remarks are drawn in Sec. IIVI 



II. THE LOCALIZATION PROCEDURE 

The new projective Wannier-Boys scheme we want to 
present here sets out from the original Wannier-Boys 
(WB) localization procedure which is discussed in detail 
in Ref. 0. Like all localization schemes for composite 
bands, it relies on the initial specification of a fixed set 
of energy bands. These bands determine the space of the 
Bloch functions which are allowed to participate in the 
multi-band Wannier transformation, the so-called active 
space. For example, the valence bands of a system can 
be chosen as such a set of bands. The WB algorithm 
is a combination of two steps: the so-called 'wannieriza- 
tion' and a Foster-Boys localization of the obtained WFs 
within the reference unit cell. Recently, the algorithm has 
been extended to operate with a multi-cell Foster-Boys 
localization to better preserve the space group symmetry 
of the system under consideration^ 

The wannierization step starts from a set of trial Wan- 
nier functions Lui n \r) which are linear combinations of 
atomic-like orbitals 



(0) 



(r) 



R) 



(1) 



We follow the notation from Ref. here. Thus, \x runs 
over all atomic basis functions <^> m (t" — Su) in the refer- 
ence unit cell, denotes their centers, and R runs over 
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all lattice vectors of the underlying Bravais lattice. To 
reduce the spacial extent of each of these WFs, the or- 
bital coefficients c^ s are set to zero for all sites s M + R at 
which the atomic Mullikcn populations 2 !^ of the given 
Wannier function falls below a certain threshold (for de- 
tails see 

Ref.0). 

The WFs obtained this way, the so- 
called 'model functions'^ are transformed to /c-space, 
projected onto the active space spanned by the selected 
BWs, orthonormalized again, transformed back to real 
space, and moved back into the reference unit cell (if 
necessary). The resulting (real) WFs, cJi , then enter 
the Foster-Boys step, where they are subject to an or- 
thogonal transformation 

N 

W «=$>«O ts (2) 

t=i 

that minimizes the spread ^[{ws 1 '}] given by the func- 
tional 

N 

0[{ W 4]=^((^|rV s }-<^|r| Ws ) 2 ) . (3) 

s=l 

Here, N is the number of energy bands involved in the lo- 
calization. Finally, the optimized functions Wj 1 ' are used 
as new (orthonormal) trial functions for the wannicriza- 
tion and the whole procedure is repeated until conver- 
gence is reached. The discarding of orbital coefficients 
and the subsequent projection onto the active space is 
the crucial part of the WB algorithm. It is combined 
with a Foster-Boys localization to ensure localization of 
the WFs also inside the unit cells. 

The described algorithm performs well for isolated 
band complexes. In the case of entangled bands, how- 
ever, the selection of proper bands to set up a suitable 
active space becomes problematic. Avoided and symme- 
try allowed crossings between the energy bands in mind 
and other, disturbing energy bands occur, and the orbital 
character one is looking for is spread over several BWs 
which in turn exhibit more or less strong admixtures from 
other, contaminating orbitals. To overcome this difficul- 
ties give up the concept of an rigid active space, abandon 
the constraint that the number of selected Bloch waves 
per fc-point has to coincide with the number N of Wan- 
nier functions per unit cell, and allow additional BWs to 
be included in the active space at each A-point. 

The selection of an appropriate A-dimensional active 
subspace is then done in the projection step during the 
wannicrization. To this end the 'model Bloch waves' 

U{r) = $> ife *6(r - R) , (4) 
ft 

which are the fc-space transformed of the model functions 
£ s (r), are projected onto the active space via 

£fc = X>«*<tfnk|U) fors = l,...iV (5) 

n=l 



with the number of selected BWs ip n k at each k- 
point being at least as large as the number of WFs per 
unit cell. The matrices U ns (k) = (ip n k\£,sk) showing up 
here establish a generalization of the unitary hybridiza- 
tion matrices used in conventional multi-band Wannier 
transformations, 2 ^. The new projected functions £f sk are 
those functions in the active space which resemble the 
initial model Bloch waves £ s fc the most. They span the 
iV-dimensional subspace used in the subsequent multi- 
band Wannier transformation. In this sense, the proce- 
dure outlined here is very similar to the one proposed by 
Souza et al}- In particular, our extended projection step 
during the wannerization can be regarded as an analog 
to the band preselection scheme used in their method. 
As for any hybrid orbitals, the orbital energies 

E-* = (Ck\F\Ck) (6) 

of the projected model Bloch waves have little in common 
with the canonical band energies e n k they originate from. 
Yet, by diagonalizing the subblock 

Fst(k) = {C k \F\e tk ) (7) 

of the Fock operator F of the system, TV new, so-called 
disentangled energy bands rj s k can be obtained. Where 
the contamination of the canonical BWs with orbitals of 
wrong character is small, the disentangled bands will es- 
sentially coincide with the canonical ones. Close to band 
crossings, where the contamination is larger, they will 
deviate substantially from the canonical bands in order 
to be able to follow the chemical nature of the underlying 
BWs. By this the disentangled bands form an effectively 
isolated complex of N bands with none of them showing 
any kinks and cones (see Sec. IIIII for more details). The 
eigenvectors resulting from the diagonalization of F st (k) 
can be regarded as sort of optimal Bloch wave hybrids 
with minimal orbital contamination and vanishing off- 
diagonal terms in the Fock operator. We will refer to 
these hybrids as disentangled Bloch waves. 

The canonical Bloch waves ip n k , to be included in the 
active space, can be selected in various ways, for exam- 
ple, by specifying an energy window and taking all BWs 
whose band energies e n k fall into this window. Alterna- 
tively, a so-called 'energy tube' around a given pair of 
reference bands (n,n) may be used, i.e. all BWs with 
band energies 

e n k 6 [enk - e, e^k + e] (8) 

are considered where e and e are some user-specified en- 
ergy tolerances. 

We have implemented the above projection and redi- 
agonalization scheme as an extension to the original WB 
localization routine in the Crystal 200x code^ a pre- 
cursor of the most recent public version of the Crys- 
tal program packaged Its ability to disentangle energy 
bands properly will be demonstrated in Sec. IIIII where 
our method is applied to the virtual bands of trans- 
polyacetylene, diamond and silicon. 
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III. RESULTS AND DISCUSSION 

All band structures shown here are calculated on the 
Hartree-Fock level of theory. The periodic ab initio pro- 
gram package Crystal (version 200x)2£l is used for that 
purpose. The localization of the WFs is performed a 
posteriori with the WB algorithmic as implemented in 
Crystal 200x 30 in conjunction with our extension for 
entangled bands which has been built into this version of 
Crystal. 

In all cases, we focus on the first few, low-lying virtual 
bands of our systems. Because of the larger extent of 
the localized virtual Wannier functions compared to the 
occupied ones, the former WFs are quite sensitive to the 
number of /c-points in the Monkhorst-Pack gridi^i We 
chosed the grids sufficiently fine to remove any ambigui- 
ties resulting from the discrete A-space integration. 



A. trans-Polyacetylene 

trans-Polyacetylene (tPA), -[HC=CH] 00 -, suits per- 
fectly as illustrative example for band disentanglement, 
because, in the basis set employed here, it exhibits three 
low-lying entangled virtual bands which are separated 
from the rest of the unoccupied band structure (Fig. ^) . 




FIG. 1: Virtual Hartree-Fock bands of trans-polyacetylene in 
the energy range from to 10 eV calculated with a cc-pVTZ 
basis set^ 

The computational parameter for tPA are taken from 
a recent study on the correlated valence and conduction 
bands of tPAtf That is, Dunning's correlation-consistent 
polarized valence triple-^ basis sets (cc-pVTZ)22 for hy- 
drogen and carbon (without / functions) are used, the 
bond distances are d(C-C) = 1.45 A, d(C=C) = 1.36 A, 
and d(C-H) = 1.09 A, and the bond angles are Z(C-C=C) 
= 121.7° and Z(C-C-H) = 118.2°. Two different geome- 
tries are considered, flat tPA (the experimental structure) 
and a distorted configuration with the hydrogens being 
bent out-of-plane by 20° in such a way that the inversion 
symmetry of the polymer is preserved. The Brillioun 
zone is sampled by a uniform grid of 100 /c-points. 

The first three virtual bands are selected for the disen- 
tanglement (Fig.GJ. One is a it* band formed by C=C it 
anti bonds, the other two are of a symmetry and describe 
C-H anti bonds (not C-C a anti bonds, as one might 



think at first glance). For the flat polymer the symme- 
try separation is perfect, for the distorted structure some 
mixture between a and it bonds occurs. Nevertheless, the 
two types of BWs remain quite different in their orbital 
character which should facilitate the band disentangle- 
ment significantly. In this sense, our first system very 
much resembles the one chosen by Souza et alw^ They 
used copper which exhibits a d band manifold which is 
entangled with a single sp valence band. 

Of course, for tPA one could localize the virtual Bloch 
waves by means of the original WB algorithm. But what 
we want to demonstrate here, is that it is also possible to 
localize the tt and a bands separately. We first consider 
the flat tPA chain. In that case, the BWs come in two dif- 
ferent symmetries and the disentanglement could simply 
be achieved by a proper labeling of the energy bands and 
the associated Bloch waves. Yet, such symmetry classi- 
fications are hard to implement in localization schemes 
for periodic systems, and thus usually not exploited. Our 
band disentanglement algorithm, however, is able to rec- 
ognize the different symmetries and to separate the bands 
properly. 

As an initial guess for the 7r-type WFs, anti-phase 
linear combinations of 2p z atomic orbitals at neighbor- 
ing carbon atoms are used. For the a-type WFs, anti- 
bonding linear combinations of 2sp 2 hybrid orbitals on 
carbon and Is atomic orbitals on hydrogen are con- 
structed. Because of this choice, the hybridization matri- 
ces U ns (k) become 3x1 and 3x2 matrices, respectively, 
with a maximum-rank subblock and all other entries be- 
ing exactly zero, as is confirmed numerically. 
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FIG. 2: Zoom into the first three virtual bands of trans- 
polyacetylene along the A line (the endpoints of the plot cor- 
respond to 1/3 and 2/3 of the TX distance) for the flat (a) 
and the distorted (b) geometry of the polymer. Solid lines 
show the canonical bands, symbols represent the disentangled 
bands: one of tt* character (•), two of a* character (o). 
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As seen in the upper panel of Fig. [21 where the canon- 
ical band energies e n k are compared to the disentangled 
band energies rj s k, our band disentanglement procedure 
is perfectly able to describe either the tt or the a orbitals 
alone. The crossing of the a bands on the left-hand side 
is an avoided one (because of the lack of any extra sym- 
metry in the small point group of the /c-points inside 
the Brillouin zone). The apparently extremely weak cou- 
pling of the two bands is due to 'soft symmetry selection 
rules'^ That means, the impact of the C-C bond length 
alternation in tPA on the C-H anti-bond dominated a* 
bands is so small that they essentially behave as if there 
would exist an additional glide plane symmetry in the 
system (like in equidistant tPA). This concept is corrob- 
orated by a detailed analysis of the involved Bloch waves. 

A more interesting situation for band disentanglement 
arises when the symmetry of the system is lowered by 
moving the hydrogen atoms out of plane. Now, the three 
energy bands avoid each other and the underlying Bloch 
waves carry both, n- and a-type atomic orbital contribu- 
tions. Using band disentanglement we are able to follow 
this contributions individually. Setting out from an ini- 
tial guess for 7r-type WFs, a single smooth band can be 
generated starting at the third canonical band at the T 
point, passing the avoided crossing in the middle of the 
Brillouin zone without any kinks and wiggles and ending 
at the lowest band at the X point (see black dots in the 
Fig. [2b) ■ The same holds for the two disentangled band 
of predominantly a character. When going from T to X 
they smoothly switch from the lower two bands to the 
upper two bands without being influenced by the com- 
plex structure of the canonical bands close to the multiple 
avoided crossing. 

After having demonstrated that our projective 
Wannier-Boys scheme is able to separate energy bands 
appropriately, the effect of the disentanglement on the 
locality of the resulting Wannier functions should be ad- 
dressed. For that purpose we turn our attention to the 
more complex case of bulk materials, diamond and silicon 
in our case. 



B. Diamond 

Because of the rather diffuse, atom-optimized basis 
functions present in the original carbon cc-pVTZ basis 
set of Dunning^ (outermost s andp exponents of 0.1285 
and 0.1209, respectively) it can not be used for a Hartree- 
Fock calculation of bulk diamond. Hence, the outermost 
exponents were reoptimized by minimizing the Hartree- 
Fock energy per unit cell of diamond. The resulting ex- 
ponents are 0.2011 for the s function and 0.6256, 0.3243 
for the p functions, typical values for diamond^ In ad- 
dition the two d functions of the triple-^ basis set were 
replaced by the single one of the corresponding double-^ 
basis set (with exponent 0.55)^ and, as for tPA, the / 
function had to be skipped because Crystal cannot han- 
dle them. This basis set, referred to as bulk-optimized 



cc-pVTZ, has been used very successfully in our embed- 
ding studies of wave-function-based correlation calcula- 
tions for diamond. 18 The experimental lattice constant of 
3.57 AiS is adopted which corresponds to an interatomic 
C-C distance of 1.546 A, together with a 40 x 40 x 40 
Monkhorst-Pack grid. 
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FIG. 3: Virtual band structure of diamond from a canon- 
ical Hartree-Fock calculation (solid lines) in comparison to 
the band energies resulting from localized Wannier functions 
(dotted lines). The localization is either done with the orig- 
inal WB algorithm (a) or by band disentanglement with an 
energy window of — 38 eV (b). 

Because the first four low-lying virtual bands of dia- 
mond are separated from the rest of the virtual bands by 
a small gap, it is possible to use the original WB algo- 
rithm to generate localized Wannier functions io s . The 
energy bands obtained from the Wannier representation 
F st (R) = (uj s (r)\F\u!t(i — R)) of the Fock operator — or 
equivalently from F st (k) as defined in Eq. — exactly 
reproduces the canonical bands (see Fig.^). Neverthe- 
less, the resulting localized Wannier functions possess 
rather substantial tails at the second nearest neighbor 
carbon atoms as seen in Fig. 0Ji. These tails can very- 
well spoil the performance of any scheme which relies on 
the loc ality of virtual WFs (like the one having been used 
in Ref.GjJ. 

Extending the active space used in the WB algo- 
rithm opens the possibility for making the Wannier func- 
tions more compact. Switching to an energy window of 
— 38 eV the spread of the virtual Wannier function, as 
measured by the Foster-Boys functional ©, reduces by 
more than a factor of two (see table P) and the undesired 
orbital contributions at the second nearest neighbor car- 
bon atoms disappear in the contour plot shown in Fig.^t. 
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FIG. 4: Contour plots of the virtual Wannier function of di- 
amond shown in the (110) plane of a C-C zig-zag chain as 
obtained by the original WB algorithm (a) or by band dis- 
entanglement with an energy window of — 38 eV (b). The 
values of the contours are ±0.046, ±0.10 and ±0.22 bohr" 3/2 
(a geometrical progression with q — 10 1//3 ). 



No shrinking of the virtual WF around the central bond 
is discernible in Fig. 0] compared to the WF obtained by 
the original WB localization. Apparently, the reduction 
in the spread is solely due to the fading of the tails in the 
localized Wannier function, precisely what we are aiming 
for. 

The price to pay, is an overall upwards shift of the dis- 
entangled bands with respect to the canonical ones, more 
pronounced at the X point than at the T or L point, 
a phenomenon we also hit on for silicon in Sec. IIII CI 
It is a feature one often observes when Foster- Boys- like 
schemes are used to localize virtual orbitals. The local- 
ization functional f2 tries to minimize the extent of the 
orbitals as much as possible, regardless of the chemical 
nature of the orbitals and, in particular, their orbital en- 
ergies. But compactness of orbitals usually implies high 
kinetic energies. Thus, it can easily happen that the WFs 
pick up more and more kinetic energy during the iterative 
procedure of the projective WB localization algorithm as 
soon as the BWs spanning active space allow so. In prac- 
tice, a compromise has to be found between tracing the 
proper orbital character in BWs energetically far away 
from the bands in mind and the risk of opening channels 
for spurious orbitals compression. 



C. Silicon 

Silicon is the last and most interesting example which 
is discussed here, because local virtual Wannier functions 
can not be generated at all for silicon without band dis- 
entanglement. 

The relativistic energy-consistent Ne-core pseudopo- 
tential from Stuttgart together with a decontracted 



TABLE I: The spread Q (per orbital) of the virtual Wannier 
functions of diamond and silicon (in bohr 2 ) as a function of 
the energy window (in eV) used for the band disentanglement. 
For comparison, the scaled spreads fi/a 2 with a being the 
lattice constant are given as well. 





diamond 






silicon 




window 


Q 


n/a 2 


window 


a, 


n/a 2 


none 


10.10 


0.792 


0-15 


33.03 


1.120 


0-38 


4.94 


0.388 


0-31 


15.45 


0.524 



[3s3p] version of the corresponding optimized valence 
double-£ basis set are used here. The basis set is aug- 
mented by a single d polarization function with an expo- 
nent of 0.4. The Si-Si distance is set to 2.352 A which cor- 
responds to a lattice constant of 5.432 A. These compu- 
tational parameters originate from the first pioneer study 
on a rigorous determination of the correlation energy of 
silicon by means of an incremental expansion, 4 ' and has 
successfully been used from that time on in all ab ini- 
tio study of correlation effects in bulk silicon performed 
with the incremental schem o 37 ! 38 or its extension to va- 
lence and conduction bands. 7 8 Like for diamond, a very 
dense 40 x 40 x 40 Monkhorst-Pack grid is used here to 
resolve the subtle details of the silicon conduction bands 
to be discussed below. 




FIG. 5: The fourth and fifth conduction band of silicon in the 
vicinity of the two symmetry-allowed crossings on the E line 
drawn as a 2-dimensional function in the rXF' plane. The 
canonical bands are shown as energy surfaces, the uppermost 
disentangled band (whose energies are only available on the 
Monkhorst-Pack grid) as dots. An energy window of — 15 eV 
is used for the disentanglement. 

The reasons why band disentanglement is absolutely 
crucial for silicon are the two symmetry-allowed cross- 
ing of a fourth and fifth conduction band on the £ line 
from X over U = K to V (the S ± S line to be precise) 
which prevents a direct application of the WB algorithm 
(see Figs. El and 0). The localization simply fails be- 
cause the active BWs exhibit symmetries different from 
the ones of the model Bloch waves £ s fc and the projec- 
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tion step JSJ) yields linear dependent projections £' sk . A 
d band is involved, one might speculate, but closer in- 
spection of the corresponding BWs reveals that the fifth 
conduction band is an sp band, formed — in contrast to 
the other four conduction bands — by s and p orbitals of 
the next atomic shell following the 3s/3p valence shell, 
at least in the basis set employed here. 

Sometimes it is argued that the band crossing problem 
discussed above can simply be solved by a proper rela- 
beling of the energy bands. This is not the case. Energy 
bands of bulk materials are three-dimensional functions 
of the crystal momentum k and the two critical bands 
exhibit an interesting topology around the symmetry al- 
lowed crossings. They only touch, twice, like the tips of 
two cones, as is clearly seen in Fig. \5\ where the band en- 
ergies are plotted as 2-dimensional energy surfaces over 
the rxr' plane. The symmetry-allowed crossings are sin- 
gularities. No band crossing occurs anywhere else in the 
neighborhood of these points. Nevertheless, the charac- 
ter of the associated Bloch waves still switches from one 
band to the other when passing from the left-hand to the 
right-hand side of the plot plane precisely as in conven- 
tional avoided crossings. 



(a) 


o (S-y'i [ -' •} 
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FIG. 6: Contour plots of the virtual Wannier function in sil- 
icon shown in the (110) plane of a Si-Si zig-zag chain as ob- 
tained from band disentanglement with an energy window of 
— 15 eV (a) or — 31 eV (b). The value of the contours are 
±0.022, ±0.046 and ±0.10 bohr" 3/2 (see Fig.gJ. 

Band disentanglement is the only way around. Two 
different energy windows are considered, one being sort 
of minimal with the upper edge at 15 eV which is closely 
above the top of the fourth conduction band, the other 
being big enough to follow the BWs with proper sp char- 
acter up to the 12-th unoccupied band at T' = — (1, 1, 1). 
In both cases, the localization could be performed with- 
out any problems. As clearly seen in Fig. [5] the different 
symmetry of the fourth and fifth conduction band along 
the E line is perfectly recognized by the band disentan- 
glement algorithm — like for the flat tPA chain discussed 
in Sec. IIII Al — and the uppermost disentangled energy 



band exclusively follows the flatter of the two canonical 
bands without any kinks. Leaving the high-symmetry 
line, the disentangled bands start to interpolate between 
the two canonical bands and form a smooth and well- 
behaved energy surface. 

The resulting localized virtual WFs of silicon are de- 
picted in Fig. Compared to diamond which exhibits 
typical sp 3 hybrid character around the nuclei the vir- 
tual WF of silicon is more symmetric and p-like in the 
vicinity of the nuclei. The same holds for the maximally 
localized virtual WF of silicon shown in Ref. [l9|. Expect- 
edly, the larger energy window yields the more compact 
WF, and like for diamond, the reduction in the spread 
is impressive, from 33 down to 15 bohr 2 (table UJ). Yet, 
all together, the virtual orbitals of silicon remain more 
diffuse than those of diamond even if the difference in the 
lattice constants is accounted for, as done by the scaled 
spreads listed in Table [I] This is not surprising, because 
of the much smaller (direct) band gap of silicon (7.3 eV 
at r, experimentally 39 ) compared to diamond (3.4 eV— ). 
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FIG. 7: Virtual Hartree-Fock band structure of silicon (solid 
lines) and the energies of the disentangled bands (dots) using 
either a — 15 eV energy window (a) or a — 31 eV energy 
window (b) for the band disentanglement. 

As already discussed above, the disentangled conduc- 
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tion bands have the tendency to shift upwards with in- 
creasing size of the energy window, a trend which is also 
found for silicon (Fig. 01. Therefore, a further increase of 
the energy window is not very helpful anymore, although 
even more compact virtual Wannier functions can be gen- 
erated this way. 

There exist a couple of further interesting features in 
the disentangled conduction bands of silicon. One is, that 
the second disentangled band along the E line climbs up 
to the upper part of the avoided crossing in the middle 
of the E panel. Apparently, close to this crossing the 
contaminating 4s/4p orbital character is solely sitting on 
the energetically more stable second Bloch wave while 
the corresponding valence 3s/3p contributions form the 
fourth BW, an interpretation which is corroborated by a 
detailed analysis of the involved BWs. 

The second point is the unexpected discontinuity in the 
uppermost disentangled band along the E line discernible 
in Fig.0L for the — 15 eV window. Its position coincides 
with the position of the avoided crossing between the 
sixth and eighth conduction band which shows up in the 
middle of the E panel between 15 and 18 eV. During the 
band disentanglement, the BWs from the upper part of 
this crossing are excluded from the active space while 
those of the lower part are present, and it seems that 
there is still a substantial part of the necessary orbital 
character present in this upper BWs to cause the abrupt 
change in the uppermost disentangled band. Obviously, 
the 3s/3p valence orbital character moves up much higher 
into the unoccupied band structure of silicon than one 
might expect at first glance. 

This is the reason why we switched to the larger en- 
ergy window of — 31 eV, though a non-negligible up- 
wards shift of the disentangled virtual bands arises. Nev- 
ertheless we consider the disentangled bands and Bloch 
orbitals from the enlarged energy window to be the more 
appropriate ones. 

IV. CONCLUSIONS 

An extension of the Wannier-Boys localization algo- 
rithm for periodic systemsi£*A£ is developed which al- 
lows to generate localized Wannier functions in the case 
of entangled energy bands. The method has been im- 
plemented into the localization routine of the Crystal 
program packaged Its main feature is the use of an en- 



larged set of active Bloch waves during the optimization 
of the unitary hybridization matrix for the multi-band 
Wannier transformation. This allows the inclusion of all 
Bloch waves which contain noticeable admixture from or- 
bitals with the same chemical character as the localized 
Wannier functions one is looking for. The proper identifi- 
cation of these admixtures is done by a simple projection 
technique. 

The efficiency of our projective Wannier-Boys algo- 
rithm is demonstrated for the virtual bands of three dif- 
ferent systems, ircms-polyacetylene, diamond and bulk 
silicon. Localized virtual Wannier functions could be gen- 
erated in all three cases. The spatial extent of them is 
found to be controllable by the size of the active space, 
i.e., the number of selected Bloch waves per fc-point. The 
more Bloch waves are considered, the more compact the 
localized Wannier functions become. Yet, at the same 
time, an increasing tendency for an overall upwards shift 
in the energies of the disentangled bands is observed. 

The same trend is discernible in the Kohn-Sham en- 
ergies of silicon discussed by Souza et al.^ though, be- 
cause of the tight energy window employed there (up to 
~11 eV), the effect is not very pronounced. The spread of 
the maximally localized virtual Wannier function of sili- 
con reported there is 30.13 bohr 2 (based on a 10 x 10 x 10 
Monkhorst-Pack grid) which is quite close to the value of 
33.03 bohr 2 we found for the — 15 eV energy window. 

The choice of the systems considered here was not ac- 
cidental. Subsequent use of the localized Wannier func- 
tions in wave-function-based post Hartree-Fock correla- 
tion methods for periodic systems which explicitly ex- 
ploit the local character of virtual Wannier functions was 
the driving force of this work, for example the method 
used in our study of the valence and conduction bands 
of irans-polyacetyleneii or in analogue investigations of 
the band structure of bulk materials like diamond^ The 
more compact the Wannier functions are the better these 
local correlation methods perform. 
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